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A density-matrix renormalization group (DMRG) method for highly anisotropic two-dimensional 
systems is presented. The method consists in applying the usual DMRG in two steps. In the first 
step, a pure one dimensional calculation along the longitudinal direction is made in order to generate 
a low energy Hamiltonian. In the second step, the anisotropic 2D lattice is obtained by coupling in 
the transverse direction the ID Hamiltonians. The method is applied to the anisotropic quantum 
spin half Heisenberg model on a square lattice. 
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The density-matrix renormalization group method in- 
troduced by White JJ a decade ago has proven a remark- 
able efficiency in the cornputation of the ground-state 

IE Q I thermodynamic [sl and dynamic 1^, |3| prop- 
erties of one-dimensional interacting electron models as 
well as electron-phonon systems [Sl llDj . An extension of 
the method to two dimensional systems is currently the 
object of some efforts But an efficient 2D algorithm 
is still lacking. 

An interesting situation which has been so far over- 
looked is the possible application of the DMRG tech- 
nique to strongly anisotropic 2D systems. In such sys- 
tems, the interaction Hamiltonian can be written as the 
sum of two terms: H = -f gH±; where TJy is a the 
sum over one-dimensional (ID) Hamiltonians (longitudi- 
nal direction) and H± is the interaction between these 
ID systems (transverse direction). H^^ and H± are of the 
same magnitude and g <C I . This type of situation arises 
in various physical systems including quasi ID organic 
and inorganic materials, weakly coupled ladder systems, 
carbon nanotube ropes, quasi ID antifcrromagnets, quasi 
ID spin-Peierls systems, etc. It is the purpose of this let- 
ter to show that a DMRG algorithm can be implemented 
for the study of such anisotropic models. 

The standard approach to anisotropic Hamiltonians 
such as H is to first diagonalize H\\ and then use its eigen- 
states to compute the corrections introduced by H± . We 
present below an algorithm which follows this spirit. This 
method consists in applying a two-step DMRG method 
to Hamiltonians having the form of H. The standard 
DMRG is first used along the longitudinal direction until 
the scale of gH± is reached. Then a transformation is 
made in order to keep only a small number ( a few tens) 
of low lying states of H^. These states are then used 
as a starting point for a second DMRG calculation along 
the transverse direction. The validity of this procedure is 
discussed. As in the thermodynamic algorithm 5] which 
was introduced by the authors, the accuracy depends on 
not only the number of states m kept by blocks, but also 
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FIG. 1: Sketch of the two-step density-matrix renormalization 
group algorithm 



on the number of the target states M. 

The two-step DMRG is illustrated in Fig. □ The first 
step of the quasi ID algorithm is the usual ID DMRG 
method. For the sake of completeness, we briefly recall 
the main steps of the algorithm. For more details, we 
refer the reader to the more extensive presentations in 
ref (3| . An iteration of the DMRG algorithm with the in- 
finite system method proceeds as follows: (a) the Hamil- 
tonian for the superblock I-I-2-1-3-I-4 ( where the blocks 
I and 4 come from previous iterations and 2 and 3 are 
new added ones. Blocks 2 and 3 are usually single site 
blocks) is diagonalized in order to obtain the few M low- 
est lying wave functions |\E'fe(I, 2, 3, 4)), fc = I,M. (b) 
the reduced density matrix of blocks I-I-2 is constructed. 
This reduced density matrix is p = ^j^cukPk, where pk 
and LUk are respectively the reduced density matrix of 
the k*^ state and its weight, pk is related to the su- 
perblock wave function \E'i;(I, 2, 3, 4) through the rela- 
tion, pfc(I, 2; I', 2') = ^3 4 1*^(1, 2, 3, 4)) (^^(I', 2', 3, 4)|. 



2 



p is diagonalized and the m eigenstates with the highest 
eigenvalues Aq. are kept. Since Aq = 1, because p is a 
density matrix, the error made by truncating this sum to 
m states is = 1 — X]a=i m ^a- ('^) These states form 
a new reduced basis in which all the operators have to 
be expanded and the block 1+2 is renamed as block 1. 
New block 2 and 3 are added to form a new superblock. 
This procedure is repeated until one reaches the desired 
lattice size. This means typically, once the energy sep- 
aration between the eigenstates of the ID system are of 
the same magnitude as g. 

Once the size L + 2 of the superblock is large enough, 
i.e when the difference between the first excited state 
and the ground state is such that \Eo — Ei\ « g, we 
perform a block transformation, of the type of the old 
RG method [HEl, on blocks 1 and 4. The two L/2- 
lattices with m states are coupled, then renormalized to 
a Z/-lattice with m states. These m states are the ground 
state and the low-lying excited states of the L-lattice. 
They can now be coupled through g. One can note that 
the boundary condition problem in the block RG method 
raised by White and Noack does not exist here, this 
is because the block transformation is only applied once. 
A 2D lattice is then generated by using the ID density- 
matrix renormalization group algorithm with chains as 
units instead of sites. 

This new block 1 will serve, in the second DMRG step, 
as a building unit for a new series of iteration that are 
similar to those performed in step 1. This block is la- 
belled block 2 in the second step since it plays the role 
played by a single site in step 1. This two-step method 
can be justified as follows: in the study of the ground 
state properties of a quasi ID system, since 5 <C 1, 
gH± may be neglected in the early steps of the renor- 
malization process. The reason is that finite size effects 
act as temperature. When the size is small, the excita- 
tion energies are far greater than the transverse coupling, 
g <^ \Eo — En\, therefore g is irrelevant. 

We now illustrate the two-step DMRG method in the 
case of weakly coupled spin-half Heisenberg chains on a 
square lattice. The Hamiltonian reads, 

Hspins — 

^ SijSj+i^( + ^ Si^iSi_i+i (1) 

i,l i,l 

where the S^^; are the usual spin-half operators. Fol- 
lowing the prescriptions above, we start by applying 
the ID DMRG on the single chain Hamiltonian TJy ; = 
J2i It was shown by White (l| that a relatively 

small number of m states kept in block 1 can lead to an 
astonishingly high accuracy for the ground-state energy 
and spin correlation functions. An accuracy ranging from 
10~^ to better than lO^^*' for the ground state energy of a 
L — 28 sites chain was obtained by keeping only to = 16 
to TO = 64 states with modest computational power. 



m 




Ladaer i 


Laader Z 


2U 


16 


-13.98792 


-13.62493 


-14.00113 


-14.07066 


24 


-13.99650 


-13.69023 


-14.01248 


-14.08803 


32 


-13.99724 


-13.70638 


-14.01526 


-14.09810 


40 


-13.99729 


-13.71446 


-14.01749 


-14.10489 



TABLE I: Ground state energies per chain of single chain 
(ID) with L = 32, of a two-leg ladder 2 x L obtained with the 
usual DMRG (Ladder 1), of the same ladder obtained with 
the two-step method (Ladder 2), of a 32 x 32 system (2D). 

A set of these highly accurate low energy eigenstates 
\<pni > and eigenvalues e„, of H\\ i is first obtained by 
the ID DMRG. The eigensets of the disconnected chain 
Hamiltonian = H\\ i are given by the relations: 

$[„] - Y[<»\^n, >, = (2) 

[n] = (ni, rt2, f^^ corresponds to an eigenset on 
the chain I. The eigenset <&[n]i diagonalizes H\\: 

< $[„]|i7|||$[„] >= £;[„]^[„],[™]. (3) 

The matrix elements of H± in the basis formed by the 
are: 

; i 

This last equation 2|can readily be transformed to: 

<$Hii?^i$H>=EsD'""s::;™^- (5) 

ii 

where S"J'™' |Si,i I'/'m, >• The renormalized 

Hamiltonian H in the basis of the <!>[„] is thus: 

^«5I^nI*W ><*Hl + ^±Es,MS.,i+i. (6) 

[n] il 

H is studied by the usual ID DMRG algorithm, ft is 
however to be noted that the computational requirement 
for H are greater than in the study of a simple spin chain 
performed in the first step. The renormalized spin oper- 
ators Si^i are now to x to matrices instead of 2 x 2. The 
size of the superblock in the second step is to'^ x m? (it 
is 4 X TO^ in step 1) if to' states are kept in the two ex- 
ternal blocks. We were nonetheless able to perform the 
calculations on a workstation for up to to = 48. 

It is instructive to compare the results which are ob- 
tained by the proposed method with a known case. A 
simple situation is the two-leg spin ladder model. We 
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FIG. 2: The binding energy of the 32 x Z system (relative to 
the isolated 32 site chain) for m = 32 and J± =0.1 




-0.446, 
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FIG. 3: The ground state energy per site for the 2D lattice 
for J± — (open circles) and J± =0.1 (filled circles) as a 
function of L for m — 32. 



m 


A(L) 


5E{L) 


dimSz ~ 


16 


0.11671 


0.58664 


12 251 


24 


0.11780 


0.65044 


66 850 


32 


0.11777 


0.72977 


205 920 


40 


0.11774 


0.78363 


477 670 



TABLE II: Finite size spin gap A(L), band- width of the 
retained states 5E{L) and the dimension of the subspace with 
= in the two-step method. 



have computed the ground state energy of the ladder 
system in tvi^o ways. First, we use the usual DMRG al- 
gorithm which have been widely applied for the study 
of this system during the last decade 0]. Second, we 
compute the ground-state energy of the two leg ladder 
by first obtaining the low energy Hamiltonian of a single 
chain. For this, we targeted Af = 5 states corresponding 
to the lowest states of the spin sectors Sz = ±2, ±1 and 
0. It is necessary to target other spin sectors when build- 
ing the low energy Hamiltonian because, when forming 
the tensor product, they may lead to lower energy states 
than some Sz = states. A criteria for the selection of 
spin sector to be selected is similar to the one applied in 
the thermodynamic algorithm. If the lowest state of Sz 
is higher than the highest state retained in Sz = 0, the 
subspace Sz is not retained. We then coupled two low 
energy ID Hamiltonians with J± to form a ladder. This 
corresponds to performing the second step of the two- 
step DMRG with one iteration and two blocks instead of 
four. 

In tabled we compare the ground state energies of the 
ladder obtained by the two methods as well as those of 
a single chain for L — 32. The infinite system DMRG 
method with open boundary conditions was used. The 
energy per site of the 32 site single chain for m = 40 is 
—0.4374 which is still far from the exact ground state en- 
ergy of the infinite chain —0.4431. This is because chains 
with open boundary conditions converge very slowly (as 



l/L) to the thermodynamic limit. It can be seen that the 
energy of the two step DMRG are systematically better 
than those of the simple DMRG. The latter are found 
to be even higher than those of a single chain. This 
can be understood by the fact that the DMRG method 
is less accurate in situations where the ground state is 
nearly degenerate. This type of situation is, for instance, 
encountered in the small coupling region of the Kondo 
lattice where a large number of states must be kept 
to get reliable results. Another look at this is that when 
J± is small, the two chains are nearly independent, the 
m states are shared between the two chains -y/m states 
for each chain . 

The ground state energy of the two dimensional sys- 
tem (32 X 32) is shown in the last column of tableland 
Fig.|21shows the variation of the ground state energy per 
chain as a function of the chain number the reference 
is a single chain with L = 32. It displays a rapid conver- 
gence as / increases. We kept the same number of states 
in the two external blocks as in the internal blocks during 
the second step of the DMRG. Truncation errors are be- 
tween 10~^ and 10~^ for m = 16 to 40. One should note, 
however, that these errors are relative to the truncated 
Hamiltonian. Hence they do not necessarily represent 
the accuracy of the method. It is important to verify 
that the truncated Hamiltonian generated for the single 
chain is accurate enough, for the ground state and for 
the low lying state, to be used as a building block for 
the 2D lattice. This will depend indeed on m and M, 
but also on J± and L. One can easily see that for a 
fixed L and J± <C A„{L), Aa-{L) is the finite size spin 
gap, the interchain matrix elements will be negligible. 
The system will behave as a collection of free chains even 
if J± is turned on. Now if J± ^ dE{L), 5E{L) is the 
width of the retained states, the matrix elements of the 
states having higher energy which have been truncated 
out have a non-negligible contribution. We thus expect 
the two-step DMRG technique to be reliable when J±^ is 
within the range A^(L) ~ < 5E{L). Values of Aa{L) 
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FIG. 4: The averaged ground state correlation function 
C'{i,j) = (-1)1^-^1 < SiSj > of a single chain (full) and 
of the 32 X 32 system (dashed) for m = 32 and J± — 0.1 



and 5E{L) are displayed in table [HJ An estimate of the 
extrapolated ground state energy per site for the 2D sys- 
tem (Fig. 13) is -.4455 for Jj_ = 0.1 when m = 32 (the 
same value of m gives -.4428 for J± =0). This is to be 
compared to the quantum Monte Carlo (QMC) estimate 
-.4485 for the infinite system obtained using clusters with 
periodic boundary conditions [T^ . 

The possible occurrence of long-range order is studied 
by observing the behavior of spin correlation functions. 
This is however complicated by the use of open bound- 
ary conditions. It is known that the OBC induce a 
strong odd-even alternation in the spin correlation func- 
tions. This alternation is also observed for the 2D system 
with J± = 0.1. In order to avoid this alternation, we de- 
fined the spin correlation function as the average between 
these two states. This is done by computing < SiS/ > 
in two ways. Si is taken as the origin of a strong link 
or as the origin on a weak link. This averaging process 
reduces the even-odd alternation as shown in Fig. 01 for 
the single chain; the alternation is suppressed in the 2D 
system. The suppression of the alternation in the 2D sys- 
tem after averaging is probably due to the fact that the 
valence bond picture is not valid for the 2D system, even 
with open boundary conditions. The observation of the 
alternation before averaging is simply reminiscent of the 
starting point which is a single chain. 

We have shown that the DMRG can successfully be 
applied to highly anisotropic 2D systems. We used the 
spin-half Heisenberg anisotropic model on a square lat- 
tice as a test case. Even in this simple case, the compu- 
tational power required increases significantly with the 



number of block states kept as shown in the last column 
of Table, m An analysis of the program reveals that 
75 percent of the CPU time is spent during the calcula- 
tion of the superblock wave function. This step involves 
mostly matrix- vector products. A significant improve- 
ment in speed is thus expected by the use of parallel 
computers. Furthermore, the initial guess for the su- 
perblock wave function was a random vector at each it- 
eration. Another improvement of the program is to use 
a transformed state from the previous iteration as done 
in the standard DMRG method 0. 

Part of this work was performed when S. M. was at 
the Croupe de Physique des Solides, Universite Paris VI. 
We wish to thank Professor Claudine Noguera and her 
group for their hospitality. We thank J.V Alvarez and A. 
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